function y = GAM_inv (x, p, p0)

   ## usage:  y = GAM_inv (x, p, p0)
   ##
   ## inv. probit transform of gamma function

   y = normcdf(x) ;
   y = (y-p0) ./ (1-p0) ;
   y = gaminv(y, p(1), p(2)) ;
#   y(!isfinite(y)) = inf ;
   y(isnan(y)) = 0 ;

endfunction
